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Abstract: We examine general physical parameterisations for viable gravitational models 
in the f(R) framework. This is related to the mass of an additional scalar field, called 
the scalaron, that is introduced by the theories. Using a simple parameterisation for the 
scalaron mass M(a) we show there is an exact correspondence between the model and 
popular parameterisations of the modified Poisson equation fj,(a, k) and the ratio of the 
Newtonian potentials r](a,k). However, by comparing the aforementioned model against 
other viable scalaron theories we highlight that the common form of fi(a, k) and 77(0, k) in 
the literature does not accurately represent f(R) behaviour. We subsequently construct 
an improved description for the scalaron mass (and therefore fi(a,k) and rj(a,k)) which 
captures their essential features and has benefits derived from a more physical origin. We 
study the scalaron's observational signatures and show the modification to the background 
Friedmann equation and CMB power spectrum to be small. We also investigate its effects 
in the linear and non linear matter power spectrum-where the signatures are evident -thus 
giving particular importance to weak lensing as a probe of these models. Using this new 
form, we demonstrate how the next generation Euclid survey will constrain these theories 
and its complementarity to current solar system tests. In the most optimistic case Euclid, 
together with a Planck prior, can constrain a fiducial scalaron mass M = 9.4 x 10~ 30 eV at 
the ~ 20% level. However, the decay rate of the scalaron mass, with fiducial value v = 1.5, 
can be constrained to ~ 3% uncertainty. 
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1. Introduction 



Precision cosmological observations [1, 2, 3, 4, 5] performed over the past decade provide 
overwhelming evidence that the Universe is currently undergoing a period of accelerated 
expansion. As a result the so-called 'standard cosmological model' has emerged, which 
attributes the acceleration to a small but finite vacuum energy. Whilst this model is both 
simple and extremely successful, the fine tuning implicit in ensuring that the vacuum energy 
only dominates at low redshift has led to a proliferation of alternative 'dark energy' models, 
which attempt to alleviate the fine tuning by introducing additional dynamical fields [6]. 

Another approach is to postulate that gravity deviates from its General Relativistic 
(GR) description on cosmological distance scales. Perhaps the most common approach to 
modifying GR involves introducing additional fields, which mediate gravity in addition to 
the standard spin-2 graviton. In this work we consider one such class of models, f(R) 
gravity, which is described by the following action 



gd 4 x 



R + f(R) 

16ttG 



(1.1) 



where f(R) is an unspecified function of the Ricci scalar, and C m is the Lagrange density 
of any non-gravitational fields. These models introduce an additional scalar field (hereafter 
referred to as the scalaron), which is gravitationally coupled to ordinary matter. This scalar 
field has the unique property that its mass depends on the 'background' energy density via 
the so-called chameleon effect, a fact that proves important in evading solar system tests 
of gravity (see for example [7, 8, 9, 10]). A significant body of literature has been devoted 
to models of the form (1.1), and we direct the reader to [11, 12, 13, 14] for reviews on the 
subject. 

In principle one can use an arbitrary function f(R) in the action (1.1), however the 
vast majority of f(R) models have been ruled out on either theoretical or observational 
grounds. Only a small subclass of functions are still considered 'viable'. In this work we 
consider the subset of models that can be written as expansions around General Relativity 
for all R > -Rvac, where i? va c = 12HqQ a is the cosmological vacuum curvature at the 
present time. Examples of viable f(R) functions are [8, 15] 



f(R) = -m 

f(R) = XR st 
both of which can be expanded as 
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for R 3> -Rvao where A and s are the modified gravity parameters nontrivially related to 
the model specific parameters m 2 ,ci,C2, R s t,n,d (see also [16, 17, 18, 19, 20, 21, 22, 23, 24, 
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25, 26]). The fact that these functions reduce to expansions around a constant in regions 
of high curvature allow them to evade solar system tests and also reproduce the standard 
background cosmology with radiation and matter dominated epochs [27, 28]. 1 

Observational constraints on f(R) models place lower bounds on the mass of the scalar 
field, at the curvature scale for which the observation is performed. In this work we will be 
concerned with cosmological constraints that can be placed on f(R) models and associated 
modified gravity parameterisations. One can broadly classify cosmological probes in two 
categories: distance based probes (such as SNIa), which utilize the background evolution 
of the Hubble parameter, and structure based probes (such as weak lensing and galaxy 
clustering) , which require knowledge of the evolution and scale dependence of density per- 
turbations. 

The paper proceeds as follows: 

In section 2 we review the background evolution and perturbation equations for f(R) 
models under the so called 'quasi-static' approximation and explore a simple parameterisa- 
tion to encapsulate them. In section 3 we show that the scalaron parameterisation that we 
initially describe is equivalent to common parameterisations of the Poisson equation fi(a, k) 
and anisotropic stress rj(a,k). Then, by comparing between scalaron models we highlight 
these forms to be less flexible than at first thought. We therefore suggest an improved 
parameterisation, which accurately represents viable f(R) models in the literature. We 
discuss potential modified gravity signals in the expansion history in 4.1, and in section 4.2 
we consider modified gravity effects in the CMB angular power spectrum and the linear 
and non linear matter power spectrum. In section 5 and 6 we examine how future CMB 
and weak lensing surveys (which probe the matter power spectrum) will constrain such 
models and argue for their complementarity to solar system constraints. We discuss our 
analysis and conclude in section 7. 

2. The Field Equations and Scalaron Mass Parameterisation 

We begin by briefly reviewing the modified gravity framework and parameterisations that 
we will be using in the paper (see [29] for further details). The f(R) field equations are 
a system of fourth order differential equations for the metric, however for models which 
reduce to expansions around General Relativity one can use the quasi-static approxima- 
tion [32, 33] to simplify them considerably. 2 The quasi-static approximation is applicable 
whenever the scalar field mass M is large relative to the background expansion ^background : 
M 2 » ^background- All current constraints on f(R) models [35, 36, 37, 38] indicate that 

lr rhe functions (1.2, 1.3) in fact still possess potential fine tuning issues in the very early Universe, see for 
example [15]. However, we only consider redshifts z < 1000 and assume that these issues are ameliorated 
by modifying the f(R) functional form at high energies. 

2 The quasi-static approximation reduces the order of the field equations, and as such we will be dis- 
carding certain solutions. However it can be shown that the additional propagating degree of freedom in a 
cosmological background is an oscillatory mode, which will decay in the early Universe and can be set to 
zero in what follows (see [15, 34] for a discussion of this point). 
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our use of this approximation scheme is justified, however we are also careful to check 
retrospectively that any parameter constraints we project for these models respect the as- 
sumption M 2 » ^background- We °-° n °t discount the possibility that more exotic models 
might be constructed which violate the quasi-static approximation, however we do not 
consider them here. 

In a cosmological context, the quasi-static approximation posits that the mass of the 
scalar field M(o) (which will be time dependent in a FRW spacetime) satisfies M 2 (a) 3> 
H 2 {a) for all a. For M 2 (a) 3> H 2 (a) the modified Friedmann equation can be written as 



H 2 ~^ + ^ + nA0+ 6i^ /R ~ 6H 2 ~ ^^r/rR " -grfR, (2-1) 
where i^GR is the standard ACDM Hubble parameter 

and in what follows, all /, /r, /rr are functions of -Rgr = 6i?GR + 12-f^GR- R subscripts 
denote differentiation with respect to R; fa = df /dR, /rr = d 2 f/dR 2 and we have defined 
/ = /+i? vac /2. Note also that we are assuming zero curvature, Q^o = 0. In the quasi-static 
approximation one can relate the scalaron mass to the functional form of / via 

M 2 = (2.3) 
ojRR 

It follows that the last four terms on the right hand side of (2.1), which constitute the effect 
of the scalaron on the background evolution, are suppressed by order 0(Hq R /M 2 ) <C 1 
relative to the standard GR contributions. In this way one expects w « — 1. 

The perturbation equations [39] are also simplified when using the quasi-static approx- 
imation. We consider the scalar perturbations of the metric in the Newtonian gauge 

ds 2 = -(1 + 2^)dt 2 + a 2 (t){l - 20) 7ij dx l dx j . (2.4) 

Considering the evolution of subhorizon modes k 3> aH, we can neglect all terms of order 
H 2 I'M 2 and take kef)/ a, kip/a » <p, tp, in which case the linearized equations reduce to 



+ (2 - 5) 

k 2 <P = -4vrG (^^p J a 2 Pm 5 m , (2.6) 
5 m + 2H5 m - 4vrG (|^fj) P™$™ = 0, (2.7) 

where K = k/aM{a). Superhorizon modes will evolve according to General Relativity, as 
they satisfy fe S uperh ^ aM{a) at all times. From equations (2.5 — 2.7) it is clear that the 
scalaron introduces a redshift and scale dependent modification to the anisotropic stress, 
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Poisson equation and subsequently the growth of over-densities. It is in this way that the 
scalaron's mass can form a natural and physical parameterisation for deviations to gravity 
in the f(R) framework. We discuss the relation to other common parameterisations in 
Section 3. 

The mass of the scalar field can take numerous functional forms, however its general 
behaviour is typically very simple. For models that reduce to expansions around a cos- 
mological constant, the mass can be written as a monotonic function that decays to the 
present. A simple functional form can be used to describe this behaviour; 



M 2 

= % (2-8) 

which contains two free parameters Mq and p. These parameters represent the value of the 
mass at the present time and the rate at which M(a) increases to the past, respectively. 
In principle one should introduce a cutoff in M(a) such that it does not diverge as a — > 0, 
however this issue is of no concern for the redshifts considered. ACDM is recovered in the 
limit Mo — > oo; where the scalaron mass becomes large and the field does not propagate. 
In this paper we write units for Mq 1 in [10 28 h _1 eV" 1 ], for comparison Hq 1 in these 
units is Hq 1 = 7500[10 28 h _1 eV -1 ]. We take the boundary of the quasi-static limit to 
be Mq 1 < STSflO^h" 1 eV _1 ], and ACDM is recovered if we set M _1 = 0. The only 
constraint that is imposed on the parameter p is that the mass M 2 must grow to the past 
at a faster rate than the background curvature; R/M 2 — > for decreasing a. This forces 
one to choose p > 2, to ensure that H 2 /M 2 decays during radiation domination. There 
is no upper bound and there is no reason to fix its value, however - the larger we take it 
the faster these models asymptote to General Relativity at early times. Typically one can 
consider reasonable fiducial values within the range 3 < p < 6. 

The expression (2.8) is chosen to approximately describe the evolution of the scalaron 
mass in a cosmological background. However, one should question how general this form is 
and how well it describes existing f(R) functions in the literature. In Figure 1 we compare 
the simple mass function (2.8) to a class of exponential models (dotted line) [16, 21, 22] and 
the power law model (1.4; solid line) over < z < 1. It is clear that (2.8) does not accurately 
mimic the f(R) models to which (1.4) belongs at low redshift, and it is this regime that we 
are primarily interested. This is for two reasons: firstly, cosmological surveys are mainly 
sensitive to low and intermediate redshifts; and secondly, for redshifts beyond this the mass 
of the scalaron generically grows very rapidly, and hence any differences in M(a) become 
irrelevant due to a highly suppressed modified gravity signal. Therefore, the exact form 
of f(R) and hence M(a) for higher redshifts is not particularly relevant; it is at this low 
redshift that we will observe an appreciable modified gravity signal. 

The reason for the discrepancy is that the simple function fails to take into account the 
late time A dominated epoch. During this period the Ricci scalar is approaching a constant 
future asymptote, and since M = M[R{a)] it follows that the mass will also asymptote to 
a constant in this regime. To incorporate this behaviour we modify (2.8) according to 
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Figure 1: A comparison of different scalaron masses as a function of rcdshift: simple phenomeno- 
logical (dashed; (2.8)), power law (solid; (1.4)) and exponential (dotted; [16, 21, 22]). They have 
been equated at z — 0. The phenomenological model, which is exactly equivalent to parameterised 
deviations in the Poission equation and anisotropic stress, as shown in the next section, is not 
particularly adept at matching the other models. 



where we have defined a* as the scale factor at matter-A equality: a* = (f^mo/rU) 1 / 3 - 
It is important to note that we still have two modified gravity parameters: Mo has the 
same interpretation as in (2.8) (the mass of the scalaron at the present, which must satisfy 
Mq -1 < 375[10 28 h _1 eV _1 ] in our approximation scheme) and v is the rate of increase 
of M(a) to the past, which we take to be v > 3/2. Note that M(a) is now dependent 
on the background cosmology via a*. The factors of four in (2.9) are included as the 
important scale is not a*, but rather the scale at which R transitions between matter and 
A domination. Since R oc ^Apt — 3Pj), where Pa = — PA and P m o = 0, the evolution of 
the Ricci scalar will be dominated by A at an earlier time than a*, hence the additional 
factor. This improved parameterization behaves like (2.8) in the matter dominated era (as 
expected), but has more appropriate behaviour at late times. It is this functional form 
that we will focus on for the majority of the paper. This parameterization is closely tied 
to the expansion (1.4) and therefore can be considered as an effective parameterization 
and match to viable f(R) models of the form (1.2, 1.3) considered in [8, 15] and also the 
exponential models of [16, 21, 22] in the redshift range z <1 (shown in Figure 1). 

3. Comparison to existing modified gravity parameterisations 

Before we continue to analyse the new more suitable scalaron form (2.9) and its cosmological 
consequences we first examine the overly simplified parameterised mass function (2.8) in 
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more detail. We find it is equivalent to an existing parameterisation in the literature 
which is widely used to represent modified gravity models, and highlights the limitations 
of that particular approach. We then discuss how our updated mass function M(a) can 
be rewritten in terms of an improved metric potential parameterisation for those familiar 
with that formalism. 

3.1 Parameterisation Equivalence 

In the quest to test deviations to gravity there have been many suggested methods and 
parameterisations. One popular procedure is to notice that in general the Poisson equation 
and metric potentials are altered [41, 42, 43, 63] and so write 



k 2 ip = -ATrGa 2 n(a, k)p5 m (3.1) 
^=r,(a,k). (3.2) 

In particular [42] first suggested that these two scale and redshift dependent functions take 
on forms to specifically reproduce the behaviour of f(R) models: 



1 + ftAffcV 

^ k) = iTApv (3 - 3) 

vM= l + A|fcV • (3 - 4) 

These have been further examined in [41, 43, 44] and have since become very popular. In 
this work we highlighted a simple functional form for the scalaron mass (2.8) to similarly 
represent f(R) models. In the scalaron framework the physical scale dependence is im- 
plicitly specified by the redshift dependent mass: K = k/[aM(a)]. We now exhibit the 
equivalence between these two approaches and obtain the relationship between the two sets 
of parameters (Ai, s) and (Mo,p) in (2.8). To do so we write (2.5) as 

V 3 + 2KV V3 + 2KV V ' 



1> 

so 

(3.6) 



3 + 2K 2 
iff \3 + 4K 2 / 

The right hand side of (3.6) is clearly r](a, k). Writing this out explicitly in terms of M(a): 



3 + 2** l + f^_l + H$R p - 2 



(3.7) 



one can now compare this to (3.4). Using the relations for scalar-tensor and f(R) theories 



4 A i_o a % 



/?1 = 3 = Af = 2_/32 Af (3 - 8) 



- 6 - 



which are given in [41, 43], for example, we can immediately write down 

^2 = ^1 fh = \, h = \- (3.9) 
Using these in (3.4) and equating with (3.7), we obtain the relationships 

s = 2p-2 (3.10) 
A, = ±. (3.11) 

As a check, one can do the same process for equations (2.6) and (3.1), which again renders 
the same relations. Due to (3.9) there are two free parameters in this /i(a, k) and i](a, k) 
scheme: s and Ai; and so the two parameterisations are identical. Perhaps at first glance 
this is not so surprising - the generalised forms (3.3) and (3.4) were invoked to match f(R) 
behaviour. However, we previously showed that this simple model does not adequately 
describe general f(R) behaviour and so neither do (3.3) and (3.4). 

3.2 Improved parameterisation 

We now consider how the late time A dominated epoch should be incorporated into the 
fi(a, k) and rj(a, k) parameterization. 

To this end we insert the mass function (2.9) into (2.5) and (2.6). In doing so we find 
that the A dominated epoch can be accounted for by considering the following parameter- 
ized forms 

, . 1 + f A 2 fc 2 /?(a)a si 

*** )= i + 14*w.U - (3 - 13) 

where we have introduced the function (3(a), defined as 

/?(«)= [a^ + 4a 3 ] - (Sl+2)/3 - (3.14) 

Again, there are still simply two free parameters: (Ai,si). The relationships between 
(Ai, si) and (Mo, v) in (2.9) are given by, 

8l = 6v-2 (3.15) 

This parameterisation is similar to the previous functional form at early times, however 
they differ at low redshifts due to /3(a). This is an important contribution, as it takes 
into account the fact that the mass of the scalar field typically asymptotes to a constant 
value for a ~ 1. Another advantage to using this form is its close link to physical f(R) 
models of the form (1.2)/(1.3). Current and novel approaches to modified gravity N-body 
simulations [49, 45] utilise an explicit form for f(R) and therefore with a large suite of these 
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simulations it will be possible to accurately reconstruct the nonlinear power spectrum for 
this more general gravitational parameterisation 3 . 

In what follows we will focus on the parameterization of the scalaron mass given by 
the expression (2.9). Whilst these two improved approaches are equivalent, we regard the 
scalaron mass to be a quantity with a clear physical interpretation. 

4. Modified Gravity Signatures 

The system of equations (2.1, 2.5 — 2.7) are sufficient to completely specify the background 
evolution and dynamics of the linear perturbations in viable f(R) models, which are subject 
to the quasi-static approximation 4 . We now illustrate the typical modified gravity signals 
resulting from these equations using the phenomenological model (2.9). 

As described previously viable f(R) models closely reproduce the ACDM expansion 
history, hence we expect only a very weak modified gravity signal in the Friedmann equation 
(2.1) and a significantly stronger signal in the perturbation equations. This statement can 
be justified using a simple intuitive argument. The scalar field mass M(a) introduces an 
additional scale into the field equations, in addition to the Hubble parameter H(a) and 
the scale of the perturbations k. Perturbations deep inside the horizon satisfy k 2 » a 2 H 2 , 
and we have chosen the mass such that M 2 3> H 2 . However the relative size of k and 
M is undetermined. One can think of M(a) as a 'modified gravity horizon'; modes that 
satisfy k 2 < a 2 M 2 will evolve as they would under the standard General Relativistic field 
equations (with order k 2 /(a 2 M 2 ) corrections), whilst modes which satisfy k 2 > a 2 M 2 will 
exhibit order unity deviations from GR. It follows that a significant modified gravity signal 
will only be present above a particular (redshift dependent) scale k mg = aM(a). This is 
clear from equations (2.1,2.5 — 2.7); corrections to the standard Friedmann equation are of 
order 0(H 2 /M 2 )), whilst the matter power spectrum exhibits large deviations from GR 
above a particular (redshift and model dependent) scale. We illustrate this in the next two 
subsections. 

The above argument is supported in [46], where constraints were imposed on viable 
f(R) models using distance based cosmological probes. There it was found that only very 
weak constraints can be placed on modified gravity parameters, even when using several 
combined data sets. 

4.1 The Expansion History 

We have argued that viable f(R) models closely mimic the standard ACDM back- 
ground expansion history. To confirm this, we take the general scalaron model (2.9) and 
compare the luminosity distance-redshift relations. 

By substituting the mass M(a) into the modified Friedmann equation (2.1), we can 
construct the luminosity distance, defined in the usual way 

3 Indeed these existing simulations have been for models resembling (1.2) and (1.3). 
4 Superhorizon modes are not accounted for in these equations, however they will simply evolve according 
to General Relativity. 
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Figure 2: Left Panel: The luminosity distances for a ACDM background expansion history (solid 
line; M -1 = 0, v = 1.5) and the scalaron model (dashed line; M^ 1 = 375[10 28 lr 1 cV -1 ], v = 
1.5; indistinguishable). Right Panel: The fractional difference between the ACDM cosmology and 
parameterised deviation. It is clear from both panels that this framework closely mimics ACDM. 



* w - (,+ "f s^r (4 - 1) 

The results are exhibited in Figure 2, taking the fiducial parameters that we use throughout 
the paper: a rn0 = 0.25, Q b0 = 0.05, h = 0.7, ln[A s ] = ln[2.34 x TO -9 ], ra s = 1 and r = 0.9. 
We contrast Mq 1 = 0.0 and v = 1.5 (ACDM) with M _1 = 375[10 28 h- 1 eV" 1 ] and v = 1.5 
(the quasi-static limit), to show the maximum difference. We observe practically no signal 
for the modified gravity parameters chosen. This should come as no surprise; we have 
argued above that there will only be a relatively weak signal in the background cosmology 
for any model in the quasi-static approximation scheme. We can therefore effectively fix 
w = — 1 and w a = for the duration of the paper, denned as w(a) = w + w a (l — a). 
As we will see in the following section however, for the same parameter choices there is a 
significant modified gravity signal in the matter power spectrum. 

4.2 The Growth of Perturbations 

We now study the effect of the scalar field on the CMB angular power spectrum and matter 
power spectrum. In order to do so we modify CAMB [64] to incorporate the generalised 
modified gravity perturbation equations (2.5 — 2.7), and use this to construct the modified 
source terms in the Boltzmann equations. 

One might first expect to observe a modified gravity signal in the angular power spec- 
trum on large scales, due to the change in the metric potentials 4> + ip that modify the ISW 
effect. However, as we show in Figure 3 the effect is not particularly significant (solid — > 
dashed line). The absence of a large modified gravity signal in the low / regime is due to 
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Multipole (I) 

Figure 3: Left Panel: The CMB angular power spectrum for ACDM (solid line; M^ 1 = 0, v = 1.5) 
and general scalaron (dashed line; Mq 1 = 375[10 28 h _1 eV _1 ], v = 1.5) parameterisation calculated 
with our modified CAMB. The modified gravity signal is not significant, even over large scales, 
because the scalaron is too massive. When the quasi-static assumption is relaxed the scalar field 
affects the larger scales and a noticeable modification to the ISW effect is seen (dotted line; Mq = 
1.5 x 10 4 [10 28 h- 1 eV" 1 ], v = 1.5). Right Panel: The fractional differences compared to the ACDM 
cosmology. 



the fact that f(R) models only modify gravity over a particular range of scales. On very 
large scales the scalaron is too heavy to cause significant deviations from GR, and similarly 
over very small scales where the density dependent suppression occurs. It is clear from our 
figures that the ISW effect occurs on scales that are too large to be affected significantly by 
the scalar field, at least within the reasonable bounds provided by the quasi-static regime. 
If we consider larger values of M^ 1 > i^o" 1 , the scalar field begins to affect the large scales 
(solid — > dotted line) . However this latter plot is for illustrative purposes only, as for these 
values our approximation scheme and hence field equations are not applicable. 

In addition to the CMB angular power spectrum, we solve (2.7) for <5 m to obtain the 
linear matter power spectrum. We normalise P(k, z) to its General Relativistic form at 
z = 100, where the modified gravity contributions are negligible. Performing these steps 
and choosing the same parameters as in section 4.1 (so as to compare the modified gravity 
signal between distance and growth based observables) , we exhibit the linear matter power 
spectra at z = in the left panel of Figure 4 (top pair) . We now observe a significant effect 
below a particular lengthscale, as expected. The modification in the power spectrum is 
also a redshift dependent effect. As z increases the power spectrum returns to its General 
Relativistic form as a result of the scalaron's mass increasing to the past. This is highlighted 
by the bottom pair of lines calculated at z = 1.5 in the same panel. 

The lack of signal in the CMB and expansion history compared to the matter power 
spectrum demonstrates the importance of structure based probes, such as weak lensing, in 
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Figure 4: Left Panel: Linear matter power spectra for ACDM (solid line; M^ 1 = 0, v = 1.5) 
and scalaron (dashed line; Mq 1 = 375[10 28 h~ 1 cV _1 ], v = 1.5) cosmologies. The modification 
to gravity causes a sizeable scale dependent effect in the growth of perturbations. The rcdshift 
dependence of the scalaron can be seen by comparing the top and bottom pairs of power spectra 
evaluated at redshifts z — 0.0 and z = 1.5, respectively. Right Panel: The environmental dependent 
chameleon mechanism can be seen in the mildly non linear regime. We exhibit the fractional 
difference (P(k) — Pgr(^)) / -Pgr(&) between the f(R) and GR power spectra for the model (2.9) 
with parameters M^ 1 = 375[10 28 h~ 1 cV _1 ] and v = 1.5. The dashed lines represent linear power 
spectra (P(k) and Pgr(A;) calculated with no higher order effects) and the solid lines are the power 
spectra calculated to second order. We see that the nonlinearities decrease the modified gravity 
signal. This is a result of the chameleon mechanism. The top set of lines correspond to z — and 
the bottom to z — 0.9; demonstrating that the modified gravity signal dramatically decreases for 
larger z. This is due to the scalaron mass being much larger at higher redshifts. Furthermore, non 
linear effects are less significant for increasing z. 



constraining /(i?)-type models. 
4.3 Non-linearities 

There is a vast quantity of useful information in the non-linear regime. However, this sector 
is one that can bias and invalidate results if not treated carefully. We now briefly discuss 
approaches to analysing nonlinear physics in the context of modified gravity. 

Broadly speaking, we can split the dynamics of the perturbations into three regimes. 
On large scales one can simply linearise the perturbation equations, where it is not necessary 
to yet account for non linearities in the fluid equations. Over an 'intermediate' range of 
scales one can consistently treat non-linear effects by calculating the density perturbations 
to higher order in a perturbative expansion. We first consider this approach, following the 
recent work [47] where this problem has been studied in detail. On the smallest of scales 
however perturbation theory is no longer applicable and one must resort to constructing 
fitting formulas for the matter power spectrum, calibrated with N-body simulations [48, 
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49, 50]. 

To understand the power spectrum at slightly higher k than the linear theory allows, 
and to see the effects of the scale dependent chameleon mechanism, we calculate P(a, k) 
to second order using the closure relations [47]. We leave the mathematical details for the 
interested reader in the Appendix. We solve the associated system of equations for the 
scalaron model for the two contrasting limits in the gravitational parameters: Mq 1 = 0.0, 
v = 1.5 and Mq 1 = 375[10 28 h -1 eV _1 ], v = 1.5. The remaining cosmological parameters 
are set to those used in the previous subsections and throughout the paper. We evolve 
the relevant equations over the range k = (0.001, 0.13)hMpc -1 , however the perturbative 
solution will approximately break down before k = 0.13 hMpc -1 . At precisely what scale 
the second order solution ceases to be accurate is not clear and will be model and redshift- 
dependent. However, we expect it to be of order ^nonlinear ~ 0.1 hMpc -1 at z = 0. 

In the right panel of Figure 4 we illustrate the modified gravity signal in the matter 
power spectrum; specifically we plot the fractional difference between the f(R) and GR 
power spectra (P(k) — Pgr(&)) / -Pgr(&) at redshift z = (top pair) and z = 0.9 (bottom 
pair). The dashed lines represent both P(k) and -Pgr(&) i n the linear regime, and the solid 
lines represent both using the second order calculation described in the Appendix. We 
observe that the nonlinear terms reduce the modified gravity signal, and have a tendency 
to return the power spectrum back to its General Relativistic form. This effect is essentially 
due to the backreaction of the density perturbations on the mass of the scalar field; the 
scalar field's mass is not just sensitive to the background energy density, but also the 
perturbations when one considers nonlinear terms in the fluid equations. This behaviour 
will not in general be captured in current linear to non linear mappings such as HALOFIT 
[65]. 

The previous Figure shows that nonlinear effects for these models are extremely small 
in the regime 0.001 hMpc -1 < k < 0.1 hMpc -1 and at redshift z = 0.9; they consti- 
tute an order < 0(2%) effect to the power spectrum with the parameter choices Mq 1 = 
375[10 28 h -1 eV -1 ], v = 1.5. In section 6 we consider potential constraints on these models 
from future weak lensing surveys, and the above calculation allows us to consistently use 
the linear matter power spectrum on scales k < 0.1 hMpc -1 . 

Although nonlinear effects are relatively small for the ranges under consideration, and 
although the perturbative method does not propagate far into k, perturbation theory is 
not without use as it can be used to calibrate or test methods used in the fully non linear 
regime. One approach to obtaining the power spectrum on scales k > O.lhMpc -1 is to use 
the Parameterized Post Friedmannian (PPF) fitting formula for f(R) models, developed 
in [51] 

p/ , n _ fno-cham(a, k) + c n i£ 2 (a, k)P GR (a, k) 
U ] ~ " ~ l + c nl £ 2 (a,fc) ^ ~ " 

where -P n o-cham is the power spectrum obtained if one were to assume that General Relativ- 
ity is modified on all scales (that is, neglecting the nonlinear effects which drive the power 
spectrum back to its General Relativistic form). c n \ and X 2 (a, k) are functions that can be 
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obtained by comparing full non linear N-body results to the semi-analytic calculation of 
P(a, k) in the mildly non linear regime. 

As shown in [47], the fitting function (4.2) is extremely accurate on scales k < lhMpc -1 , 
but only if we have knowledge of -P n o-cham m the nonlinear regime from simulations. Whilst 
promising work has been made in this direction [45, 48, 49, 50], there are currently insuf- 
ficient N-body results to calibrate the fitting function for different scalaron parameters. It 
is for this reason we do not use (4.2) in the next section of the paper, as it does not guar- 
antee a better fit to f(R) models than halofit in the absence of an accurate description 
of -Pno-cham- Moreover, an industrious programme of N-body simulations to explore the 
gravitational parameter space is unlikely to be embarked upon unless a suitable, general 
and physical parameterisation/framework is settled upon. Addressing this is one of the 
intentions of this work. 



5. Weak Lensing 

Weak lensing is a statistical measure of the distortion of distant source galaxies. It occurs 
due to gravitational shear induced by the intervening mass distribution, and so observing 
the weak lensing signal yields information on the distribution of matter projected along 
the line of sight. 

With respect to analysing cosmological models, of particular interest is the convergence 
k. This is related to the shear matrix $^ via the relation 

$11 + $22 

k= (5.1) 

where the shear matrix is in turn related to the Newtonian potential 4> as 

rxh 

*ij = / dxg(x)W (5-2) 
J o 

and 

g( X ) = 2 [ Xh dA(xO r(x)r K~ X) - (5-3) 

where r{\) is the comoving angular diameter distance and n(x) is the probability of finding 
a galaxy at comoving distance \- The resulting convergence power spectrum is given by 
[66] 

where P(k, \) is the matter power spectrum and Xh is the comoving distance to the horizon. 

It is clear therefore that weak lensing is sensitive to the underlying matter power 
spectrum and the expansion history. This makes the probe highly capable of distinguishing 
between modified gravity and dark energy models. In addition, a significant advantage to 
using weak lensing is that the signal is sourced by the metric potentials and is therefore 
not subject to any unknown galaxy bias. Modified gravity models can also act to alter the 
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relation between the power spectrum of the metric potentials P^p and the matter power 
spectrum P$ [52, 53], which is expressed finally in (5.4). 

Weak lensing does face technical challenges in the form of stringent requirements for 
shape estimation [54] and intrinsic alignments [55] but the potential power of this method 
for dark energy, neutrino masses and testing gravity provides vast motivation for solving 
them. We focus on this probe of the matter power spectrum in the next section and 
combine it with the CMB after. 



6. Future Surveys and their Gravitational Limits 
6.1 Weak Lensing: Euclid 

We perform a Fisher matrix analysis to calculate expected parameter sensitivity for a future 
weak lensing survey: the proposed space-based Euclid mission [56] (see [40, 41] for related 
works). While one can undertake a more studious MCMC analysis for this purpose, as in 
[41], the Fisher matrix method produces bounds that are essentially equivalent for these 
illustrative purposes [41]. Euclid intends to measure galaxies over a sky coverage of 20'000 
square degrees with an effective galaxy density of 40 gal/arcmin 2 . The corresponding Fisher 
matrix is given by 

^zS^-f («■« 

where C is the weak lensing observable 

Ci j (l) = P ij + (lL) 6 f. (6.2) 

Pij is given in (5.4) with indices denoting tomographic redshift bins, and fi{ is the average 
galaxy number per steradian in bin i. The quantity (j^) 1 ^ 2 is the rms intrinsic shear 
in each component, which we set to 0.22. dC/dp\ is the derivative of the weak lensing 
observable with respect to the parameter p being varied. Cov is the covariance matrix 



Cov 



c&W'Ci*] = fwnt a/ fc&(0cji(0 + qg(/)c*J . (6.3) 



(21 + l)/ sky A/ 

and / s ky is the sky fraction. To model the redshift distribution we use the expression 

n(z) oc z a exp[-(z/z f] (6.4) 

with zq = Zm/1.412, a = 2, j3 = 1.5. z m is the median redshift and for a Euclid-like 
survey z m = 0.9 [56, 58]. We consider five tomographic redshift bins, chosen such that 
there are approximately an equal number of galaxies in each. Finally, to account for the 
photometric error we convolve the redshift bins with a Gaussian of width a z = a p (l + z), 
where a p = 0.03. 
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6.2 The CMB: Planck 



In addition to the weak lensing survey, we use a forecast for the upcoming Planck surveyor 
[59]. Similarly the CMB Fisher matrix is given by 

where Cx,i is the harmonic power spectrum for the temperature-temperature (X = TT), 
temperature-E-polarization (X = TE) and the E-polarization-E-polarization (X = EE) 
power spectra. The covariance matrix COVxy for the power spectra is given, under a 
Gaussian assumption, by, 



COV T ,T 


= fi 


{Ctj + W^B' 2 ) 2 


(6.6) 


cov E , E 


= fi 


(Cej + W^B- 2 ) 2 


(6.7) 


COWte,te 


= fi 


Ctej + (C T:l + W^B- 2 ) (C E:l + W P l B- 2 ) 


(6.8) 


COV t ,e 


— ftC TEl 


(6.9) 


COV j 1 TE 


= ftC TE ,i {C T jl + W^B~ 2 ) 


(6.10) 


CO\e,te 


= fiC TE ,i {C E ,i + W P l Bl 2 ) , 


(6.11) 



where fe = ( 2 £+i)/ k an( ^ Wt,p = (o"T,p#fwhm)~ 2 is the weight per solid angle for tempera- 
ture and polarization, with a la sensitivity per pixel of <jt,p with a beam width of #f w h m - 
The beam window function is expressed in terms of the beam full width half maximum as 
5, = ex P (-^ + l)^ hm /i61n2). 

In the forecast we act conservatively and use only one of Planck's channels (143 GHz) 
for the signal assuming the other frequencies are used in foreground removal. The channel 
has a beam of 0f w h m = 7.1' and sensitivities of cjt = 2.2/j.K/K and up = 4.2 fiK/K. We 
take the sky fraction to be / s k y = 0.80 due to Galactic coverage. 

Over the lowest multipoles one needs an accurate treatment of the Integrated Sachs- 
Wolfe (ISW) effect via perturbations. Due to uncertainty in this treatment and potential 
polarisation foreground contamination we exclude the lowest multipoles from our analysis. 
Instead we utilise the range 30 < i < 2000. A consequence of this is slightly weaker limits 
owing to a poorer determination of the optical depth to reionisation r. However, we feel 
that conservativeness is preferable to an incorrect or overly optimistic treatment. This 
cut in scale also has limited affect on the direct determination of the scalaron parameters. 
This is because, as we demonstrated earlier, the CMB is largely insensitive to gravitational 
modifications of this kind. This was shown in Figure 3. Therefore, in this analysis the 
role of the Planck fisher matrix will be to act as a prior on each of the main cosmological 
parameters and aid the limits by breaking parameter degeneracies. 

For the inclusion of the Planck Fisher matrix we start from the fiducial parameter 
set ft m h 2 , 8s, ln(.As), ilbh 2 , n s and r, where 6s is the angular size of the sound horizon 
at last scattering and ln(As) is the logarithm of the primordial amplitude of scalar per- 
turbations. This is a more Gaussian parameter set as detailed by the Dark Energy Task 
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Force (DETF)[60]. We calculate the Planck CMB Fisher matrix using camb for our power 
spectra. The optical depth acts as a nuisance parameter for the matter power spectrum 
and so we marginalise over it in our analysis. We then transform the resulting parameters 
to Euclid's basis set: U m , h, In (-As), flf, and n s . The resulting Fisher matrix can then be 
added to the weak lensing matrix. The transformation is performed using the Jacobian 




(6.12) 



and the final basis of the Fisher matrix is determined by 

F = JFJ T . (6.13) 

6.3 Forecast Results and Analysis 

For the fiducial model, we take the same parameters used throughout the paper: £l m o = 
0.25, n b0 = 0.05, h = 0.7, ln(A s ) = ln(2.34 x 10~ 9 ) and ra s = 1. The equation of state is 
effectively fixed to w = — 1 and w a = given our earlier arguments in Section 4.1 regarding 
the expansion history. We make the simplifying assumptions of zero spatial curvature 
^fco = an d that there are no massive neutrinos. We note that including massive neutrinos 
could act to introduce an additional degeneracy, as they suppress the growth of structure 
[61]. 

The ACDM limit of the scalaron parameterisation in (2.9) corresponds to M 1 — > 0. 
We wish to use this cosmology as one of the fiducial models, however if we take Mq 1 = the 
resulting observables are completely insensitive to v due to there being no modified gravity 
signal. In this case the fisher matrix is non-invertible and we will be unable to derive forecast 
constraints on the cosmological parameters. It follows that when we take ACDM as the 
fiducial model, we must temporarily fix the time evolution: v = 1.5. Another procedure is 
to choose Mq" 1 / and vary both modified gravity parameters. We therefore also undertake 
analyses with fiducial cosmologies in a modified framework: MJ" 1 = 75[10 28 h _1 eV _1 ] and 
v = 1.5. 

We perform our analyses in two regimes: in the first instance we only consider weak 
lensing contributions from linear scales 10 < I < 400. The resulting constraints that we 
obtain will be conservative but accurate, as we have a firm understanding of the pertur- 
bations in this range. For the second set of analyses we include the full range of scales 
10 < I < 10000, however in doing so we include non linear effects. In this work we use 
the HALOFIT fitting formula to reconstruct the non linear power spectrum. We argued 
earlier that with the present limited number of simulations this was at least as good as 
the new formalism given by (4.2). We note that this approach will not yield an accurate 
description of the nonlinear regime for these modified gravity models at present, however 
the resulting forecasts will be indicative of the gain that could be obtained if we had a 
suitable description of the power spectrum on small scales. 

The limits on these models are illustrated in Figure 5 for Euclid only and Euclid and 
Planck combined. We have also included the constraints arising from solar system tests of 
these models. A discussion of the solar system bound is included in the Appendix. 
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Figure 5: Top Left Panel: Euclid-only projected bounds on the scalaron with fiducial values 
Mq 1 = (ACDM) and v = 1.5 (fixed). The wider set of limits arc for linear scales (£ max = 400) 
and the tighter set are with all scales included (£ ma x = 10,000). In all plots solid (dashed) lines 
represent 68% (95%) confidence limits. Top Right Panel: Similarly the improvement from using 
non linear data is evident when both parameters are varied (Mq 1 = 75[10 28 h _1 eV _1 ], v = 1.5 
fiducial): £ max = 400 (less constraining) — > £ max = 10,000 (most constraining). It is not possible 
to currently treat this regime accurately but the improvement is indicative of the gain possible. 
Both sets of contours are for Euclid + Planck. Bottom Left Panel: The specific improvement from 
including Planck (tighter set) can be seen when added to Euclid-only full non linear limits (wider 
set). Bottom Right Panel: The combined limit from the previous panel is shown in contrast to 
current solar system constraints (dotted line; region above excluded). The models are successful 
in evading such tests (note the log vertical axis) and so weak lensing surveys will provide essential 
and far superior bounds. 



It is clear that weak lensing will be most adept at probing the nature of gravity: The 
Euclid probe is capable of stringent constraints on the general gravitational parameterisa- 
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tion introduced and examined in this work. However, Figure 5 highlights the significance 
in understanding the non linear regime in this gravitational framework. The additional 
information, particularly where the scale dependent signal is evident, leads to vastly re- 
duced bounds and could be vital in making a definite measurement of new gravity. This is 
exhibited clearly in the top right panel. 

The bottom left panel illustrates the great importance of complementarity in cosmol- 
ogy. Despite the fact that Euclid is essential in providing the actual modified gravity signal, 
Planck is highly beneficial in reducing the bounds. This contribution is an indirect one 
and it occurs by reducing the parameter space. Specifically, Planck tells us the primordial 
amplitude A s and spectral index n s of the ~ unmodified power spectrum (Figure 3) with 
exquisite detail. Therefore any subsequent change to the shape or, approximately, the am- 
plitude of the matter power spectrum, as in Figure 4, must be the result of alterations in 
gravity through [i and v. Indeed, this is exactly what we see; Planck breaks the degener- 
acy mainly through the spectral index and, to a slightly lesser extent, the amplitude. This 
explains the much tighter combined limit in the bottom left panel of Figure 5. 

The overall degeneracy between the modified gravity parameters for lensing only can 
simply be understood by considering that as the scalaron's mass decreases there is a larger 
modified gravity signal and a boost to the power spectrum. To compensate this and to give 
the same observable one can increase the time evolution through v such that the scalaron 
decays back to GR faster in the recent past. Finally, the bottom right panel in Figure 5 
demonstrates the absolute necessity in using cosmological probes to explore gravity and this 
general parametersation. The density dependent suppression makes it severely challenging 
even for precision solar system tests. 

7. Conclusions 

We have examined a general scheme for a class of viable gravitational f(R) models; the 
scalaron framework. This is physically motivated and gives rise to distinctive observational 
features. For viable f(R) models the background expansion is nearly identical to ACDM, a 
fact reflected in a luminosity distance-redshift relationship that is practically indistinguish- 
able from the standard cosmology. By solving the modified perturbation equations with 
an augmented CAMB, we also find that the CMB does not add information on the modified 
signal directly, because the scalaron mass is too large to affect the ISW substantially. The 
most significant observational effect for these models occurs in the matter power spectrum; 
specifically P(a, k) acquires a redshift and scale dependent modification relative to GR. 
This signal can be used to constrain modified gravity parameters with weak lensing. 

By examining parameterisations of f{R) gravity in detail we have argued that existing 
functional forms in the literature do not capture the essential features of f(R) models at 
late times z < 1. To accurately describe the behaviour of these models for z < 1 we have 
directly parameterized the scalaron mass M(a). Whilst our approach can be related to the 
standard fi(a,k) and r/(a,k) parameterization in the literature [41, 42, 43, 63] (and also 
the 7 parameterization; see [29, 30, 31]), the mass of the scalar field is an ideal function 
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to parameterize, in the sense that it is independent of scale, only weakly dependent on 
cosmological parameters and has a clear physical interpretation. 

Using the general scalaron model (2.9), we found that a future weak lensing probe 
such as Euclid will be able to place extremely tight constraints on deviations to gravity. 
We quote our final results as M = 1.34 ± 0.62 x l(T 30 [heV], v = 1.50 ± 0.18 for I < 
400 and M = 1.34 ± 0.25 x 10~ 30 [heV], v = 1.50 ± 0.04 for the I < 10000 (Planck + 
Euclid) combined analysis. Our forecast results highlight a number of interesting points, 
specifically; the importance of constraining these models with cosmological weak lensing 
data (compared to solar system tests), the value of using CMB data to break parameter 
degeneracies, and also the importance of modeling the nonlinear power spectrum. As 
stated previously, this regime is particularly hard to treat in any gravitational theory, and 
could bias the results if implemented blindly. While some approaches, such as higher order 
perturbation theory, can give us an indication of physics at nonlinear scales, it is evident 
that there is a huge demand for programmes of simulations and a better approach to fitting 
to them. A detailed study of this problem will be presented elsewhere. 
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8. Appendix I: Nonlinear regime 

The mildly nonlinear regime can be treated semi-analytically using a standard, order by 
order expansion of the perturbation equations. Specifically, for f(R) models one can write 
the fluid equations as [47] 

d$ a n „ NJR f d 3 k lC l 3 k 2 



dr 



(2vr)3 fe(k " kl " k 2)^c(ki, k 2 , r)* 6 (ki , r)<f c (k 2 , r) (8.1) 



/d k d k d k 
-5 D {k - ki - k 2 - k 3 )a a6cci (ki, k 2 , k 3 , r)$ 6 (ki, r)$ c (k 2 , r)$ d (k 3 , r) 

where we have used the linearized Einstein equations to remove the metric potentials, and 
have defined 
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'k 2 M 2 

H(a,k) =(- + _) (8.9) 



M 2 (a) = 3^ (8.10) 

M 3 = f? (8.11) 

and r = log [a]. By using the standard expansions <E> a = + <I>i 2 ' ) + <E>i 3 ' ) + ... one can 
solve (8.1) and obtain the power spectrum to any order using the definition 

(2vr) 3 5 D (k + kOP afe (|k|,r) = ($ a (k,r),d» 6 (k / ,r)) (8.12) 

However, rather than using this approach we instead follow [47] and directly calculate 
the matter power spectrum to second order, using the so-called closure relations. Specifi- 
cally one defines two additional statistical quantities in addition to P a &(a, k); 

(27r) 3 fe(k + k / )i?a6(|k|,r,r / ) = (^a(k,r),$ 6 (k / ,r / )) r > r' (8.13) 
\3x f u v\n nu^u/ i$ »M\ 



(2vr) 3 fe(k - k')G ab (|k|,r,r') = ^^ ) ) r > r' (8.14) 
By solving the following linearized equations for R ab (\k\,T,T r ) , G a &(|k|, r, r') 

A ab (k,r)R bc (\k\,r,r') = (8.15) 

A a6 (k,r)G fec (|k|,r,r / ) = (8.16) 

where 

A„ 6 = 5 a6 — + O a& (k,r) (8.17) 

using initial conditions i? a &(|k|, n, Tj) = P a &(|k|, Tj), G a {,(|k|, Tj, Tj) = <5 a &, one can construct 
the next to leading order corrections to the power spectrum by solving the equation 

T abcd (k,r)P cd (\k\,T) = f T dr"M as (k,T,r")R bs (k,r,T") + f N as (k,r,T")G bs (k,T, T ") 

+S as (k, r)P sb (k, t) + (a o b) (8.18) 

where 
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and 



d 

Fabcd = + Sac^bdik, t) + S bd ^ac{k, t). (8.19) 

It is the solution to the above set of equations that is used in section 4.3. 
9. Appendix II: Solar system constraint 

Finally, we discuss the origin of the solar system constraint imposed in Figure 5. Whilst 
many different constraints have been imposed on f(R) models, we feel that the solar system 
bound considered in (for example) [8] is both robust and model independent, and it is for 
these reasons that we incorporate it. 

We direct the reader to [8] for the details of the derivation, and simply state the result, 

I/rCRs^SxICT 11 (9.1) 

where R s = 87rGp s , and p s is the typical background density at which solar system tests are 
performed at. We take this conservatively to be R s = 10 5 Hq [10, 67], although a thorough 
treatment would require a detailed description of the solar system density distribution. 

It is straightforward to convert this bound into a constraint on the f(R) model pa- 
rameters for the function (1.4); we find 

/ R \ 2s+l 

\f R (R s )\=2Xs(^\ ~ 2Asl(r 5 ( 2s+1 > < 5 x 10" 11 . (9.2) 

However for the model (2.9) we have not explicitly written a functional form f{R). There- 
fore to apply the bound (9.1) we must first calculate the relationship between M(a) and 
Jr. Using the a <C 1 limit of (2.9), rearranging the expression (2.3) and integrating, we 
find 

-f on 2 f -Gu—Aj- ^mOA* 2 6u-3 o\ 

f R = — 3S2 m0 /U / a u da = -- -a u °. (9.3) 

JO Av — i. 

In addition, we require the value of the scale factor at which i? cos = R s ~ 10 5 i? vac , 

Rfi „K ^ 
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Hence a s = (O m o/4^A) 1//3 10 5//3 . Substituting this into (9.3) we find that the solar system 
constraint for the model (2.9), 

( ^mo\ ^ 1 10 -5(2„-l) < 5 x 10 -11 (9 5) 

2v -\ V4^a/ 



- 24 - 



